Distribution of relative velocities in turbulent aerosols 
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We compute the distribution of relative velocities for a one-dimensional model of heavy parti- 
cles suspended in a turbulent flow, quantifying the caustic contribution to the moments of relative 
velocities. The same principles determine the corresponding caustic contribution in d spatial dimen- 
sions. The distribution of relative velocities A.v at small separations R acquires the universal form 



p(Av,R) ~ R" 



for large (but not too large) values of |Av|. Here D2 is the phase- 



space correlation dimension. Our conclusions are in excellent agreement with numerical simulations 
of particles suspended in a randomly mixing flow in two dimensions, and in quantitative agreement 
with published data on direct numerical simulations of particles in turbulent flows. 
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Collisions of particles in randomly mixing or turbulent 
flows ('turbulent aerosols') have been studied intensively 
for several decades. The main goal, not yet achieved, is 
to find a reliable parameterisation of the collision ker- 
nel. This is an important question since the frequency 
of collisions between suspended particles determines the 
stability of turbulent aerosols. Direct numerical simu- 
lations of particles in turbulent flows (see for example 
P, HI) show that collision velocities (and thus the colli- 
sion kernel) increase precipitously as the 'Stokes number' 
St is varied beyond a threshold. This dimensionless pa- 
rameter, St = (7r) _1 , is defined in terms of the particle 
damping rate 7 and the relevant correlation time r of 
the flow - both explained below. In Q (see also [|J) this 
steep increase was attributed to the fact that singulari- 
ties (so-called 'caustics') in the particle dynamics result 
in large relative velocities at small separations. Caustics 
occur when phase-space manifolds describing the depen- 
dence of particle velocity upon particle position fold over. 
In the fold region, the velocity field at a given point in 
space becomes multi-valued, giving rise to finite veloc- 
ity differences between particles at the same position. In 
Q it was argued that the collision rate exhibits an acti- 
vated St-behaviour, exp(— A/St), reflecting the frequency 
at which caustics occur d, y-Q . 

The collision kernel is determined by the distribution 
of relative velocities of the suspended particles at small 
spatial separations. In order to calculate this collision 
kernel from the microscopic equations of motion, it is 
necessary to understand how caustics determine the dis- 
tribution of collision velocities, and to compute how caus- 
tics affect the moments of this distribution. This is an 
important question for many reasons, not only because 
the St-dependence of the collision rate is of fundamental 
importance in aerosol physics. First, such an approach 
could show which aspects of the collision dynamics in tur- 
bulent aerosols are universal, and which aspects depend 
on the particular nature of the turbulent velocity fluctu- 
ations. Under which circumstances are highly simplifcd 
models appropriate (such as a white-noise velocity field 
with a single spatial scale)? Second, it is known that 
turbulent aerosol particles form fractal clusters (see [f| 
and references cited therein). These fractal objects fold 



in phase space when the above-mentioned singularities 
occur. This raises the question: how do fractal spatial 
fractal clustering and caustic singularities together deter- 
mine relative velocities in turbulent aerosols? 

In this paper we compute the joint steady-state distri- 
bution of relative velocities and particle separations, first 
in a one-dimensional model Q for relative velocities of 
inertial particles using perturbation theory. Our results 
show that the distribution of relative velocities at small 
separations is a power law, and that the power is deter- 
mined by the correlation dimension of the steady-state 
phase-space fractal. We show that this power law is a 
consequence of caustics, and that caustics make a sub- 
stantial contribution to relative velocities at small sepa- 
rations, consistent with the picture outlined above, and 
the findings of Ref. [3j . Relative velocities in higher spa- 
tial dimensions are determined by the same principles, 
and we show that the corresponding power law (Eq. ([7]) 
below) is a universal property of turbulent aerosols. It 
determines, in general, how caustics and fractal cluster- 
ing give rise to fractal phase-space distributions which in 
turn determine the collision kernel. Last but not least, 
our findings explain the results of a recent scaling anal- 
ysis of the particle- velocity structure functions obtained 
by direct numerical simulations of turbulent flows 1 9] . 

We first analyse the one-dimensional model p| 



x = j(u(x, t) — v) 



(1) 



Dots denote time derivatives, v = x is the particle ve- 
locity, 7 the rate at which the inertial motion is damped 
relative to the fluid, and u(x, t) is a Gaussian random 
field describing the velocity of the fluid. Fig. Q] shows 
particle trajectories according to (QJ and illustrates how 
a smooth manifold of initial conditions develops fold sin- 
gularities (caustics). We assume that u(x,t) is charac- 
terised by spatial and temporal correlation scales r\ and 
r, and by its typical size, uq (the statistical properties 
of u(x,t) are described in the caption to Fig. [TJ. The 
dynamics is then determined by two dimensionless pa- 
rameters, by the Stokes number St and by the 'Kubo 
number' Ku = upr/rj which characterises the fluctua- 
tions of u(x, t) (flEI- I n turbulent flows the Kubo num- 
ber is of order unity, but in the following we analyse the 
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FIG. 1: a Particle trajectories x(t) according to Eq. ([T]) 
with a Gaussian random field u(x,t), with (u) — and 
(u(x,i)u(0, 0)) = «o ex P [— 2: 2 /(2t7 2 ) — \t\/r] for small values 
of x, and periodic boundary conditions in x with period L. 
Parameters: t — 0.01, 77 = 0.1, L = 1, «o = 1, 7 = 4/9, that 
is Ku = 0.1, St = 225 (particles initially at rest). Panels b to 
d show how the corresponding phase-space manifold folds. 

dynamics in the limit of small values of Ku where the 
suspended particles experience u(x, t) as a white-noise 
signal. In the limit Ku — > 0, the joint density pi Av, Ax) 
is determined by a Kramers equation [ill \v\ with a 
Ax-dependent diffusion constant. In dimensionless units 
it' = jt, x' = x/r/, v' = vj (777), v! — uj (777) and dropping 
the primes to simplify the notation) we have: 

dtp = -d Ax (Avp) + d Av {Avp)+V{Ax)d Av p, 

1 f°° 

V{Ax) = - di(Aw(Ax,t)Aw(Ax,0)) . (2) 

2 j-00 

Here Au(Ax,t) — u(x + Ax,t) — u(x,t). For the model 
described in Fig.Q]we have V(Ax) ~ e 2 Ax 2 for | Ax| <C 1 
and V(Ax) ~ 2e 2 = V for \Ax\ > 1. Here e 2 = Ku 2 St. 

In the limit of Ax — > 0, the steady-state solution 
of ([2|) is found by separation of the variables Ax and 
z = Av / Ax. Inserting the ansatz g^{Ax)Z^{z) with sep- 
aration constant /1 into (|2|) results in g M (Ax) = jAxI^ , 
while Zn{z) solves = — pzZ ll {z)+d z {z+Z 2 +e 2 d z )Z ll {z). 

The steady-state solution of is given by a weighted 
sum of g^Zfj, over the allowed values of \x. We know that 
the distribution of spatial separations exhibits a power 
law as Ax — > 0, corresponding to spatial clustering. We 
expect that the dominant contribution to this distribu- 
tion derives from the smallest positive allowed value of 
p.. For certain values of /u, the functions Z^z) are known 
in closed form (for p = | and for p = — 1 13]). For 
other values of /x, we determine an analytical solution by 
expanding Z^ around /1 = 0. We find: 

w = E(4) / d^...dw(ni» 

fc=0 i=li=0 
2fe+l i 

x exp(-^(-l)V(^t J )) (3) 

i=0 j=0 




FIG. 2: a Steady-state relative-velocity distribution 
p(Av,Ax) at Ai = 6 x 10~ 4 (o), 5 x 10" 2 (0), and 5 x 10" 1 
(□), otherwise same parameters as in Fig. 1. Also shown are 
results from ((3| subject to the constraint = (solid 
red lines), b Typical trajectory in the Aa;-At)-plane. Ar- 
rows indicate the direction of motion, c Results for fi c {^ 2 ) 
from ((3| , red line, numerical results for the spatial correlation 
dimension di (o), and numerical results for the phase-space 
correlation dimension D2 (□), both for Ku = 0.1. Also shown 
are results of a linearised Langevin model (Ku = 0) for cfo (V) 
and D2 (0). Two critical values of e 2 are shown: for the path- 
coalescence transition [j§], e 2 x ~ 1.77, and the second critical 
value e 2 2 ~ 20.7 where 1I2 equals unity, d Numerical results 
(symbols) for the distribution p(Vr, R) of relative velocities 
from computer simulations of a linearised Langevin model for 
inertial particles suspended in a two-dimensional incompress- 
ible Gaussian random flow, e 2 = 0.04, at R = 5 x 10 -6 (o), 
5 x 10~ 4 (0), and 5 x 10~ 2 (□). Also shown is the theoretical 
result Eq. (8) (solid red line), the value of /i c = D2 was deter- 
mined by numerically computing the phase-space correlation 
dimension D2- 

with V{y) = er 2 (y 2 /2 + y 3 /3) and t = z. A related 
expansion was used in [14[ to calculate moments of the 
finite-time Lyapunov exponent for particles accelerated 
in a random time-dependent potential. To determine the 
allowed values of fi we consider the large-z asymptotes of 
Z^(z). The known solution Z (z) || exhibits power-law 
tails ~ \z\~ 2 reflecting the fact that caustic singularities 
are reached in finite time. Correspondingly, when z is 
large we can neglect the term d z (zZ fi ) in the equation for 
Z^. The resulting equation is solved by combinations of 
Kummer functions with the asymptotic behaviour Z^ ~ 
|z| M ~ 2 for large values of \z\. A power-law ansatz in the 
equation for Z^ gives 

^(z)~A±(±^- 2 (l + ( M -l) 2 - 1 + ...) (4) 

for z — > ±00. Eq. (0| suggests that p(Av, Ax) approaches 
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power-law form as Aa; — > + 

P (Av, Ax) = \Ax\^- 2 Z^(Av/Ax) ~ A* (±Awf ~ 2 (5) 

for positive and negative values of Aw, respectively. Now 
we use that interchanging the particles in a particle pair 
(Aa; — > — Ax and Av —> —Av) cannot alter the distribu- 
tion. In particular this must be true at Aa; = which 
requires = A~ . The allowed values of fx can be com- 
puted by evaluating ^ and imposing that = A~ . To 
simplify the problem we only consider the smallest pos- 
itive allowed value of p (referred to as p c , and shown in 
Fig. Ht). This yields an accurate description (except in 
the very far tails) of the distribution p(Av, Ax) for small 
values of Aa; as Fig. shows. 

Three regions can be distinguished in p(Av,Ax). 
First, the body of the distribution (labeled 1 in Fig. [5^) 
corresponds to small values of z (\z\ -C 1 that is \Av\ <C 
|Aa;|). Typically, the separation remains constant in this 
region, Aa; Axo, and Av obeys dAv/dt — Av + Au. 
In the limit Ku — > 0, St — > oo so that e 2 = Ku 2 St remains 
constant, this is an Ornstein-Uhlenbeck process. It gives 
rise to a broad Gaussian distribution of relative veloc- 
ities, approximately independent of Au. This smooth 
contribution is not universal. 

Second, the regime \z\ S> 1 (that is |Aw| 3> |Ax|, la- 
beled 2 in Fig. [5^.) corresponds to particles approach- 
ing each other on different caustic branches. This sin- 
gular contribution, due to caustics, yields large values of 
| Aw | as | Ax | —> 0. In this case, the deterministic part 
of Eq. © dominates, giving rise to linear trajectories 
Av = Avo + Axo — Ax. When \Av\ ^> |Aa;|, we have 
Au Avq. In this limit, the distribution p(Av, Ax) 
becomes approximately independent of Aa;. As shown 
above, it exhibits the power-law form ([5]). The power-law 
is clearly visible in Fig. This singular caustic contri- 
bution to p(Av, Ax) gives rise to large moments of Aw, 
and thus to a large collision rate. We emphasise that the 
dynamics in this regime is universal, that is independent 
of the particular statistics of the velocity field. Fig. ^Bp 
shows a typical phase-space trajectory passing through 
the regimes 1 and 2 and contributing to the smooth and 
singular parts of the distribution p(Aw, Ax). 

Third, as Fig. [5h shows, the power law is cut off at 
very large velocities, ensuring that the moments of Av 
do not diverge. This region (labeled 3 in Fig. is not 
described by keeping the smallest positive value of p c 
only. The cut off, Av c , is simply due to the fact that 
the relative velocities Awo at large separations are of or- 
der T>y 2 . This suggests that Av c ~ 2?c/ 2 - The form 
of p(Aw,Ax) at | Aw | ^> Av c depends upon the details 
of the model. In turbulent flows at large Stokes and 
Reynolds numbers, Aw c is determined by the 'variable- 
range projection principle' suggested in [12j. The cut off 
Av c corresponds to a cut off Au c /Aa; in the distribution 
of z. Thus p(Ax, z) does not strictly factorise, except 
in the limit Aa; — > 0. The probability distribution of z 
for large values of \z\ is obtained by integrating the joint 
distribution over Aa; to ±Aw c /|z|. As a consequence of 



the power law ((4]), the tails of Z (z) found in [1] are re- 
covered, namely Zq(z) ~ M~ 2 - This shows once more 
that the distribution of the finite differences Aw and Aa; 
is governed by the same principles as the fluctuations of 
dv/dx. 

What do these results imply for the moments 
TO p (Aa;) = JdAv\Av\ p p(Av,Ax) at small Aa;? We find: 

m p (Ax) = IAxI^" 1 fdz\z\ p Z^{z) 

~ {Ax^+^ibp + c p |Aa;|- p ^ +1 ) . (6) 

Here the coefficient b p results from the contribution of 
the body of Z Pa (z) to the z-integral in ([5]). The sec- 
ond term is the caustic contribution, it results from in- 
tegrating the power-law tails of Z flc (z) up to the cut- 
off Aw c /Aa;, keeping only the leading-order behaviour in 
(H]). This form of m p (Ax) was deduced from results of 
numerical simulations for the distribution of relative ve- 
locities and separations in a one-dimensional Kraichnan 
model |15|, and it was suggested that the exponent p c 
equals the spatial correlation dimension. This is correct 
for /i c < 1, as the following argument shows. The spatial 
correlation dimension c?2 describes the power-law scaling 
for p — 0, namely mo(Aa;) ~ |Aa;| d2_1 as Aa; — > 0. From 
Eq. © we deduce mo(Aa;) ~ |Ax| roui ^ ,1 ^ _1 , implying 
that p c — d,2 for < p c < I . But for p c > 1 this is not the 
case. A change of variables Aa; — > Aw = V Aa; 2 + Aw 2 in 
([S]) allows us to compute the distribution of 'phase-space 
separations' Auj. We find that it exhibits a power-law of 
the form Awj Mc_1 as Aw — > 0. This proves that p c is in 
fact equal to the phase-space correlation dimension D2 
(which ranges between and 2). It also demonstrates 
that D2 = d,2 for < D2 < 1 (consistent with results 
of numerical simulations of particles in two-dimensional 
random flows at finite Kubo numbers [Hj]). We computed 
spatial and phase-space correlation dimensions from nu- 
merical simulations of the equation of motion ([TJ. The 
results are shown in Fig. [2t and are in good agreement 
with the analytical theory. Caustic contributions dom- 
inate for p > 1 — (j, c , the smooth contribution is thus 
irrelevant for p = 1. 

For two- and three-dimensional incompressible flows, 
the distribution of relative velocites A.v at spatial sep- 
aration R has a form similar to ([5]): as R — > 0, the d- 
dimensional analogue of equation ([2]) is solved by the 
ansatz g ^(R) Z ^(z) where z = Av/R. We find that 
Zfj, exhibits tails of the form Z^(z) ~ \z\^~ 2d . Com- 
parison with the expected form of the distribution of 
(phase-)space separations demonstrates that p c equals 
the phase-space correlation dimension. We find: 

p{Av,R)~ R d - x \Av\ D ^- 2d . (7) 

Now consider the distribution of relative radial velocities 
Vr = Av ■ R (R is the radial unit vector). For |Vr| ^> R 
this distribution is found by integration of ([7]). Assuming 
D2 < d+1, we find for small values of R (and |Vr| ^> R) 

p{V R ,R)~R d - 1 \V R \ D *- d - 1 . (8) 
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Eq. (|8| yields the following expression for the moments 
of m p (R) of Vr\ 

m p (R) = b p (St)RP +D2 - 1 + c p (St)i? d - 1 , (9) 

of the same form as ([6]). But note that here the singular 
caustic term acquires a geometrical factor Setting 
p = in @ shows that c? 2 = min{Z? 2 j^} in analogy 
with the one-dimensional case. For p = 1, Eq. © is 
equivalent to the ansatz proposed in [3j, apart from the 
difference that in Eq. (7) of Q, the asymptotic behaviour 
Z?2(St) — > d as St — > was used for the first term in ([9]). 

Fig. [2ji shows numerical results for a linearised 
Langevin model for inertial particles in a two-dimensional 
incompressible Gaussian random flow. Shown is the dis- 
tribution p(Vr,R). As in the one-dimensional case, the 
power-laws in ([7]) and (HJ represent the universal contri- 
bution of caustics to the distribution of relative veloci- 
ties at small separations. This demonstrates how caus- 
tics and fractal clustering determine the distribution of 
relative velocities at small separations. Our numerical 
results imply that the prefactor c p (St) of the singular 
caustic term in © is of activated form c p (St) ~ Q~ A p/ f2 
for small e (not shown). 



In the remainder we show how our theory explains the 
results of direct numerical simulations of particles sus- 
pended in turbulent flows @. In particle- velocity 
structure functions S P {R) = ni p (R)/mo(R) were com- 
puted numerically and analysed in terms of the scaling 
ansatz S P (R) ~ R^ p . Fig. 3 in [9j shows the scaling expo- 
nent £i as a function of St. Our theory, Eq. ©, predicts 
that £i ~ 1 for small values of St, and £i = d — efe for 
large values of St. These predictions are in quantitative 
agreement with the results presented in Fig. 3 in @ . Our 
theory also explains the distribution of zr = vr/R shown 
in Fig. 5 in Eq. © yields P(z R \R) ~ |2jj| i?2_d ~ 1 

which explains the power law in the inset of Fig. 5 in 
0. The far tails in this figure are expected to be de- 
termined by the variable-range projection principle [l2j 
yielding P(zr\R) ~ exp(— C \zr i?| 4 / 3 ). As a next step 
it is necessary to determine whether the above expres- 
sion describes the functional form of the far tails of this 
distribution in (l7l |. 
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